<!DOCTYPE HTML PUBLIC "-//W3C//DTD HTML 4.01 Transitional//EN"
                "http://www.w3.org/TR/REC-html40/loose.dtd">
<html>
<head>
  <title>Description of solve_am</title>
  <meta name="keywords" content="solve_am">
  <meta name="description" content="">
  <meta http-equiv="Content-Type" content="text/html; charset=iso-8859-1">
  <meta name="generator" content="m2html v1.5 &copy; 2003-2005 Guillaume Flandin">
  <meta name="robots" content="index, follow">
  <link type="text/css" rel="stylesheet" href="../../../m2html.css">
</head>
<body>
<a name="_top"></a>
<div><a href="../../../index.html">Home</a> &gt;  <a href="../../index.html">CO4OI</a> &gt; <a href="#">Solvers</a> &gt; <a href="index.html">AM_solvers</a> &gt; solve_am.m</div>

<!--<table width="100%"><tr><td align="left"><a href="../../../index.html"><img alt="<" border="0" src="../../../left.png">&nbsp;Master index</a></td>
<td align="right"><a href="index.html">Index for CO4OI/Solvers/AM_solvers&nbsp;<img alt=">" border="0" src="../../../right.png"></a></td></tr></table>-->

<h1>solve_am
</h1>

<h2><a name="_name"></a>PURPOSE <a href="#_top"><img alt="^" border="0" src="../../../up.png"></a></h2>
<div class="box"><strong></strong></div>

<h2><a name="_synopsis"></a>SYNOPSIS <a href="#_top"><img alt="^" border="0" src="../../../up.png"></a></h2>
<div class="box"><strong>function [sol, of, crit_AM]=solve_am(ymeas, param) </strong></div>

<h2><a name="_description"></a>DESCRIPTION <a href="#_top"><img alt="^" border="0" src="../../../up.png"></a></h2>
<div class="fragment"><pre class="comment"> 

 solve_AM - Solves the following minimization problem      

       min_{x,y,z} ||A(x \otimes y \otimes z)||_2^2 s.t. sol \in R_+,

 using Alternate Minimization algorithm, as a generalized version of algorithm presented in [1] 
       
   INPUT
       - ymeas: vector of measurements
       - param: data structure with parameters of the algorithm

   OUTPUT
       - sol: local solution
       - of: objective function
       - crit_AM: stoping criteria used

   param is a Matlab structure containing:

   General parameters:
 
   - verbose: 0 no log, 1 print main steps, 2 print all steps.

   - n,m: dimensions of the sought image

   - T: table indicating the probed triplets of frequencies

   - max_iter: max. nb. of iterations (default: 200).

   - rel_obj: minimum relative change of the solution (default:
   1e-4)
       The algorithm stops if
           | ||x(t)||_2 - ||x(t-1)||_2 | / ||x(t)||_2 &lt; rel_obj and
           | ||y(t)||_2 - ||y(t-1)||_2 | / ||y(t)||_2 &lt; rel_obj and
           | ||z(t)||_2 - ||z(t-1)||_2 | / ||z(t)||_2 &lt; rel_obj, 
       where x(t), y(t), z(t) are the partial estimates of the solution at iteration t.

   Parameters for the backtracking procedure in Forward-Backward Algorithm

   - L: initial approximation Lipsitsz constant

   - eta: tuning parameter (to find the optimal step size in the gradient step)

   - max_iter_BT: max. nb. of iterations (backtracking)


   [1] J.P. Haldar and D. Hernando, &quot;Rank-Constrained Solutions to Linear Matrix
   Equations Using PowerFactorization&quot;, 2009, IEEE Signal Processing Letters, 16, 584</pre></div>

<!-- crossreference -->
<h2><a name="_cross"></a>CROSS-REFERENCE INFORMATION <a href="#_top"><img alt="^" border="0" src="../../../up.png"></a></h2>
This function calls:
<ul style="list-style-image:url(../../../matlabicon.gif)">
<li><a href="solve_fb_bt.html" class="code" title="function [sol, crit_ncFB] = solve_fb_bt(y, param)">solve_fb_bt</a>	solve_fb_bt - solve the following problem</li></ul>
This function is called by:
<ul style="list-style-image:url(../../../matlabicon.gif)">
</ul>
<!-- crossreference -->



<h2><a name="_source"></a>SOURCE CODE <a href="#_top"><img alt="^" border="0" src="../../../up.png"></a></h2>
<div class="fragment"><pre>0001 <a name="_sub0" href="#_subfunctions" class="code">function [sol, of, crit_AM]=solve_am(ymeas, param)</a>
0002 <span class="comment">%</span>
0003 <span class="comment">%</span>
0004 <span class="comment">% solve_AM - Solves the following minimization problem</span>
0005 <span class="comment">%</span>
0006 <span class="comment">%       min_{x,y,z} ||A(x \otimes y \otimes z)||_2^2 s.t. sol \in R_+,</span>
0007 <span class="comment">%</span>
0008 <span class="comment">% using Alternate Minimization algorithm, as a generalized version of algorithm presented in [1]</span>
0009 <span class="comment">%</span>
0010 <span class="comment">%   INPUT</span>
0011 <span class="comment">%       - ymeas: vector of measurements</span>
0012 <span class="comment">%       - param: data structure with parameters of the algorithm</span>
0013 <span class="comment">%</span>
0014 <span class="comment">%   OUTPUT</span>
0015 <span class="comment">%       - sol: local solution</span>
0016 <span class="comment">%       - of: objective function</span>
0017 <span class="comment">%       - crit_AM: stoping criteria used</span>
0018 <span class="comment">%</span>
0019 <span class="comment">%   param is a Matlab structure containing:</span>
0020 <span class="comment">%</span>
0021 <span class="comment">%   General parameters:</span>
0022 <span class="comment">%</span>
0023 <span class="comment">%   - verbose: 0 no log, 1 print main steps, 2 print all steps.</span>
0024 <span class="comment">%</span>
0025 <span class="comment">%   - n,m: dimensions of the sought image</span>
0026 <span class="comment">%</span>
0027 <span class="comment">%   - T: table indicating the probed triplets of frequencies</span>
0028 <span class="comment">%</span>
0029 <span class="comment">%   - max_iter: max. nb. of iterations (default: 200).</span>
0030 <span class="comment">%</span>
0031 <span class="comment">%   - rel_obj: minimum relative change of the solution (default:</span>
0032 <span class="comment">%   1e-4)</span>
0033 <span class="comment">%       The algorithm stops if</span>
0034 <span class="comment">%           | ||x(t)||_2 - ||x(t-1)||_2 | / ||x(t)||_2 &lt; rel_obj and</span>
0035 <span class="comment">%           | ||y(t)||_2 - ||y(t-1)||_2 | / ||y(t)||_2 &lt; rel_obj and</span>
0036 <span class="comment">%           | ||z(t)||_2 - ||z(t-1)||_2 | / ||z(t)||_2 &lt; rel_obj,</span>
0037 <span class="comment">%       where x(t), y(t), z(t) are the partial estimates of the solution at iteration t.</span>
0038 <span class="comment">%</span>
0039 <span class="comment">%   Parameters for the backtracking procedure in Forward-Backward Algorithm</span>
0040 <span class="comment">%</span>
0041 <span class="comment">%   - L: initial approximation Lipsitsz constant</span>
0042 <span class="comment">%</span>
0043 <span class="comment">%   - eta: tuning parameter (to find the optimal step size in the gradient step)</span>
0044 <span class="comment">%</span>
0045 <span class="comment">%   - max_iter_BT: max. nb. of iterations (backtracking)</span>
0046 <span class="comment">%</span>
0047 <span class="comment">%</span>
0048 <span class="comment">%   [1] J.P. Haldar and D. Hernando, &quot;Rank-Constrained Solutions to Linear Matrix</span>
0049 <span class="comment">%   Equations Using PowerFactorization&quot;, 2009, IEEE Signal Processing Letters, 16, 584</span>
0050 <span class="comment">%</span>
0051 
0052 
0053 
0054 <span class="keyword">if</span> ~isfield(param, <span class="string">'verbose'</span>), param.verbose = 1; <span class="keyword">end</span>
0055 <span class="keyword">if</span> ~isfield(param, <span class="string">'rel_obj'</span>), param.rel_obj = 1e-3; <span class="keyword">end</span>
0056 <span class="keyword">if</span> ~isfield(param, <span class="string">'max_iter'</span>), param.max_iter = 200; <span class="keyword">end</span>
0057 
0058 
0059 n=param.n;
0060 m=param.m;
0061 
0062 iter=0;prev_x=0;prev_y=0;prev_z=0;
0063 
0064 <span class="keyword">while</span> 1
0065     
0066 <span class="comment">%% define Ayz operators to solve for x</span>
0067 <span class="comment">% compute Fourier trnasforms y and z</span>
0068 yhat=1/n*fft2(reshape(param.yinit,n,m));
0069 zhat=1/n*fft2(reshape(param.zinit,n,m));
0070 
0071 param.A=@(x)Ayz(x,yhat(:),zhat(:),param.T);
0072 param.At=@(x)Ayz_t(x,yhat(:),zhat(:),param.T);
0073 
0074 param.solini=param.xinit(:);
0075 
0076 
0077 <span class="comment">% solve problem for x</span>
0078 
0079 [x, crit_ncFBx]=<a href="solve_fb_bt.html" class="code" title="function [sol, crit_ncFB] = solve_fb_bt(y, param)">solve_fb_bt</a>(ymeas,param);
0080 
0081 <span class="comment">% update initial solution</span>
0082 param.xinit=x(:);
0083 
0084 
0085 <span class="comment">%% define Axz operators to solve for y</span>
0086 <span class="comment">% compute Fourier trnasforms x and z</span>
0087 xhat=1/n*fft2(reshape(param.xinit,n,m));
0088 zhat=1/n*fft2(reshape(param.zinit,n,m));
0089 param.A=@(x)Axz(x,xhat(:),zhat(:),param.T);
0090 param.At=@(x)Axz_t(x,xhat(:),zhat(:),param.T);
0091 
0092 
0093 param.solini=param.yinit(:);
0094 
0095 <span class="comment">% solve problem for y</span>
0096 
0097 [y,crit_ncFBy]=<a href="solve_fb_bt.html" class="code" title="function [sol, crit_ncFB] = solve_fb_bt(y, param)">solve_fb_bt</a>(ymeas,param);
0098 
0099 <span class="comment">% update initial solution</span>
0100 param.yinit=y(:);
0101 
0102 
0103 <span class="comment">%% define Axy operators to solve for z</span>
0104 <span class="comment">% compute Fourier trnasforms x and y</span>
0105 xhat=1/n*fft2(reshape(param.xinit,n,m));
0106 yhat=1/n*fft2(reshape(param.yinit,n,m));
0107 param.A=@(x)Axy(x,xhat(:),yhat(:),param.T);
0108 param.At=@(x)Axy_t(x,xhat(:),yhat(:),param.T);
0109 
0110 param.solini=param.zinit(:);
0111 
0112 <span class="comment">%solve problem for z</span>
0113 [z, crit_ncFBz]=<a href="solve_fb_bt.html" class="code" title="function [sol, crit_ncFB] = solve_fb_bt(y, param)">solve_fb_bt</a>(ymeas,param);
0114 
0115 <span class="comment">% update initial solution</span>
0116 param.zinit=z(:);
0117 
0118 
0119 
0120 <span class="comment">%% Global stopping criterion</span>
0121     
0122        
0123 rel_varx = norm(x - prev_x)/norm(prev_x);
0124 rel_vary = norm(y - prev_y)/norm(prev_y);
0125 rel_varz = norm(z - prev_z)/norm(prev_z);
0126 
0127 <span class="comment">% compute objective</span>
0128 rel_error= norm(A(x,y,z,param.T)-ymeas)/norm(ymeas);
0129 of=.5*norm(A(x,y,z,param.T)-ymeas)^2;
0130 
0131 <span class="comment">% final solution</span>
0132 sol=1/3*(x+y+z);
0133 
0134 <span class="keyword">if</span> param.verbose &gt;= 1
0135     fprintf(<span class="string">'  ||Ax-b||_2 = %e, rel_var = %e\n'</span>, <span class="keyword">...</span>
0136         norm(A(x,y,z,param.T)-ymeas), rel_varx);
0137 <span class="keyword">end</span>
0138 
0139 <span class="keyword">if</span> (rel_varx &lt; param.rel_obj &amp;&amp; rel_vary &lt; param.rel_obj &amp;&amp; rel_varz &lt; param.rel_obj || norm(x-prev_x)==0)
0140     crit_AM = <span class="string">'REL_NORM'</span>;
0141     <span class="keyword">break</span>;
0142 <span class="keyword">elseif</span> (rel_error &lt; 1e-3)
0143     crit_AM = <span class="string">'REL_ERROR'</span>;
0144     <span class="keyword">break</span>;
0145 <span class="keyword">elseif</span> iter &gt;= param.max_iter
0146     crit_AM = <span class="string">'MAX_IT'</span>;
0147     <span class="keyword">break</span>;
0148     
0149 <span class="keyword">elseif</span> strcmp(crit_ncFBx,<span class="string">'MAX_IT'</span>) || strcmp(crit_ncFBy,<span class="string">'MAX_IT'</span>) ||strcmp(crit_ncFBz,<span class="string">'MAX_IT'</span>)
0150     crit_AM = <span class="string">'MAX_IT'</span>;
0151 <span class="keyword">end</span>
0152     
0153     
0154     
0155 <span class="comment">% Update variables</span>
0156 iter = iter + 1;
0157 prev_x = x;
0158 prev_y = y;
0159 prev_z = z;
0160 
0161 <span class="keyword">end</span>
0162 
0163 <span class="keyword">end</span>
0164</pre></div>
<hr><address>Generated on Thu 18-Jul-2013 19:25:15 by <strong><a href="http://www.artefact.tk/software/matlab/m2html/" title="Matlab Documentation in HTML">m2html</a></strong> &copy; 2005</address>
</body>
</html>